Nucleation in scale-free networks 
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We have studied nucleation dynamics ol the Ising model in scale-free networks with degree distri- 
bution P(k) ~ k ' by using forward flux sampling method, focusing on how the network topology 
would influence the nucleation rate and pathway. For homogeneous nucleation, the new phase 
clusters grow from those nodes with smaller degree, while the cluster sizes follow a power-law distri- 
bution. Interestingly, we find that the nucleation rate Ruom decays exponentially with the network 
size N, and accordingly the critical nucleus size increases linearly with N, implying that homoge- 
neous nucleation is not relevant in the thermodynamic limit. These observations are robust to the 
change of 7 and also present in random networks. In addition, we have also studied the dynamics 
of heterogeneous nucleation, wherein w impurities are initially added, either to randomly selected 
nodes or to targeted ones with largest degrees. We find that targeted impurities can enhance the 
nucleation rate Rnet much more sharply than random ones. Moreover, \n(R,Het/ 'Rnom) scales as 
^7-2/t-i anc j w £ Qr targeted and random impurities, respectively. A simple mean field analysis is 
also present to qualitatively illustrate above simulation results. 

PACS numbers: 89.75.Hc, 64.60.Q-, 05.50.+q 



I. INTRODUCTION 



Complex networks describe not only the pattern dis- 
covered ubiquitously in real world, but also provide a 
unified theoretical framework to understand the inher- 
ent complexity in nature [H-Q- Many real networks, as 
diverse as ranging from social networks to biological net- 
works to communication networks, have been found to 
be scale- free i.e., their degree distributions follow a 
power- law, P(k) ~ k r . A central topic in this field has 
been how the network topology would influence the dy- 
namics taking place on it. Very recently, critical phenom- 
ena in scale-free networks (SFNs) have attracted consid- 
erable research interest M , such as order-disorder transi- 
tions 0-E3| i Percolation [Tll - ilij , epidemic spreading [l5l , 
synchronization self-organized criticality [l8l[l9| . 

nonequilibrium pattern formation (20| , and so on. These 
studies have revealed that network heterogeneity, charac- 
terized by diverse degree distributions, makes the critical 
behaviors of SFNs quite different from those on regular 
lattices. However, most previous studies focused on eval- 
uating the onset of phase transition in different network 
models. There is little attention paid to the dynam- 
ics/kinetics of phase transition, such as nucleation and 
phase separation in complex networks. 

Nucleation is a fluctuation-driven process that initi- 
ates the decay of a metastable state into a more stable 
one 21|. A first-order phase transition usually involves 
the nucleation and growth of a new phase. Many impor- 
tant phenomena in nature, including crystallization [22J, 
glass formation (23[, and protein folding [24 |. etc., are 
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associated with nucleation. Despite its apparent impor- 
tance, many aspects of nucleation process are still un- 
clear and deserve more investigations. The Ising model, 
which is a paradigm for many phenomena in statistical 
physics, has been widely used to study the nucleation 
process. Despite its simplicity, the Ising model has made 
important contributions to the understanding of nucle- 
ation phenomena in equilibrium systems and is likely to 
yield important insights also for nonequilibrium systems. 
In two-dimensional lattices, for instance, shear can en- 
hance the nucleation rate and at an intermediate shear 
rate the nucleate rate peaks |25(, a single impurity may 
considerably enhance the nucleation rate [26|, and the 
existence of a pore may lead to two-stage nucleation and 
the overall nucleation rate can reach a maximum level at 
an intermediate pore size. Nucleation pathway of Ising 
model in three-dimensional lattice has also been studied 
using transition path sampling approach [27]. In addi- 
tion, Ising model has been frequently used to test the va- 
lidity of classical nucleation theory (CNT) [HIlD]. How- 
ever, all these studies are limited to regular lattices in 
Euclidean space. Since many real systems can be mod- 
eled by complex networks, it is thus natural to ask how 
the topology of a networked system would influence the 
nucleation process of Ising model. To the best of our 
knowledge, this has never been studied in the literature 
so far. 

Although the main motivation of the present study is 
to address a fundamental problem in statistical physics, 
it may also be of practical interest. For example, our 
study may help understand how public opinion or belief 
changes in the social context [331 ] . where binary spins 
can represent two opposite opinions and the concept of 
physical temperature corresponds to a measure of noise 
due to imperfect information or uncertainty on the part 



2 



of the agent. Another example is functional transition in 
real biological networks, such as the transition between 
different dynamical attractors in neural networks [3^ |. 
wherein the two states of the Ising model may correspond 
to a neuron being fired or not. Other examples include 
gene regulatory networks, wherein genes can be on or off, 
corresponding to the two states of Ising model 35 1 . 

In the present work, we have studied nucleation pro- 
cess of the Ising model in SFNs. Since nucleation is an 
activated process that occurs extremely slow, brute-force 
simulation is prohibitively expensive. To overcome this 
difficulty, we adopt a recently developed forward flux 
sampling (FFS) method to obtain the rate and path- 
way for nucleation [36| . For homogeneous nucleation, we 
find that the nucleation begins with nodes with smaller 
degree, while nodes with larger degree are more stable. 
We show that the nucleation rate decays exponentially 
with the network size N, and accordingly the critical nu- 
cleus size increases linearly with N, implying that ho- 
mogeneous nucleation can only occur in finite-size net- 
works. Comparing the results of networks with different 
7 and those of random networks, we conclude that net- 
work heterogeneity is unfavorable to nucleation. In addi- 
tion, we have also investigated heterogeneous nucleation 
by adding impurities into the networks. It is found that 
the dependence of the nucleation rate on the number of 
random impurities is significantly different from the case 
of targeted impurities. These simulation results may be 
qualitatively understood in a mean-field manner. 

The rest of the paper is organized as follows. In SecJTTl 
we give the details of our simulation model and the FFS 
method applied to this system. In Sec lIII| we present 
the results for the nucleation rate and pathway. We then 
show, via both simulation and analysis, that the system- 
size effect of the nucleation rate and heterogeneous nu- 
cleation. At last, discussion and main conclusions are 
addressed in SeclIVI 



II. MODEL AND SIMULATION DETAILS 



A. The networked Ising model 



probability min(l, e~ l3AE ), where (3 = l/(ksT) with fee 
the Boltzmann constant and T the temperature, AE is 
the energy change due to the flipping process. Gener- 
ally, with the increment of T, the system will undergo 
a second-order phase transition at the critical tempera- 
ture T c from an ordered state to a disordered one. To 
study nucleation, we set J = 1, h > 0, T < T c , and start 
from a metastable state in which Sj = — 1 for most of 
the spins. The system will stay in that state for a signif- 
icantly long time before undergoing a nucleation transi- 
tion to the thermodynamic stable state with most spins 
pointing up. We are interested in dynamics of this nu- 
cleation process. 



B. Forward flux sampling 

FFS method has been used to calculate rate constants, 
transition paths and stationary probability distributions 
for rare events in equilibrium and nonequilibrium systems 
[lHHlllifl 38 40]. This method uses a series of interfaces 
in phase space between the initial and final states to force 
the system from the initial state A to the final state B in 
a ratchet-like manner. An order parameter \{x) is first 
defined, where x represents the phase space coordinates, 
such that the system is in state A if X(x) < Ao and state B 
if A(x) > A m , while a series of nonintersecting interfaces 
Ai (0 < i < m) lie between states A and B, such that 
any path from A to B must cross each interface without 
reaching Ai+i before A^. The transition rate R from A to 
B is calculated as 

R = $ A , P A m Ao =$ A , TT P A^+i Ai , (2) 

where <&a,o is the average flux of trajectories crossing Ao 
in the direction of B. P (X m \X ) = n£Lo P (K+i W) 
is the probability that a trajectory crossing Ao in the 
direction of B will eventually reach B before returning 
to A, and P (Ai+i|Ai) is the probability that a trajectory 
which reaches Ai, having come from A, will reach A;+i 
before returning to A. For more information about FFS, 
please turn to Ref . [4fj| . 



The Ising model in a network comprised of N nodes is 
described by the Hamiltonian 

H = - J 22 "-.r s < s .. ''X!'"" w 

i 

where spin variable Si at node i takes either +1 (up) 
or —I (down). J is the coupling constant and h is the 
external magnetic field. The element of the adjacency 
matrix of the network takes = I if nodes i and j are 
connected and otherwise. 

Our simulation is performed by Metropolis spin-flip dy- 
namics (37j . in which we attempt to flip each spin once, 
on average, during each Monte Carlo (MC) cycle. In 
each attempt, a randomly chosen spin is flipped with the 



III. RESULTS 

A. Homogeneous Nucleation: Rate and pathway 

To begin, we first consider homogeneous nucleation in 
Barabasi- Albert scale-free network (BA-SFN) whose de- 
gree distribution follows a power-law P(k) ~ A; -7 with 
the scaling exponent 7 = 3 [5j. We define the order 
parameter A as the total number of up spins in the net- 
works. We set N — 1000, the average degree (k) = 6, 
T = 2.59, h = 0.7, A = 130 and A m = 880, where T 
is lower than the critical temperature T c ~ 10.36. The 
spacing between interfaces is fixed at 3 up spins, but the 
computed results do not depend on this spacing. During 



FIG. 1: (Color online) Snapshots of nucleation in a BA scale- 
free network with N = 100 and (k) = 2 at four different 
stages. Up-spins and down-spins are indicated by red squares 
and black circles, respectively. 
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the FFS sampling, we perform 1000 trials at each inter- 
face, from which at least 100 configurations are stored in 
order to investigate the statistical properties of the en- 
semble of nucleation pathway. We obtain $,4,0 = 1-24 x 
10' 4 MC 'step' 1 spin- 1 and P (A m |A ) = 4.48 x lO" 46 , 
resulting in Rnom = 5.55 x 10~ 50 MC step' 1 spin' 1 fol- 
lowing EqjS] Such a nucleation rate is very low such that 
a brute-force simulation would be very expensive. 

From the stored configurations at each interface, one 
can figure out the details of the nucleation pathway. Fig- 
ure Q] illustrates schematically four stages of a typical 
nucleation pathway. Clearly, the new phase (indicated 
by squares) starts from nodes with smaller degree, while 
nodes with larger degree are more stable. This picture 
is reasonable because nodes with larger degrees need to 
overcome more interfacial energies. Figure [U[a) plots 
the average degree of network nodes in the new phase, 
(knew), as a function of the order parameter A. As ex- 
pected, (k new ) increases monotonously with A. On the 
other hand, it is observed that the formation of large 
clusters of new phase is accompanied with the growth 
and coalescence of small clusters. Interestingly, we find 
that the size N c of new phase clusters follows a power law 
distribution at early stages of nucleation, P(N C ) ~ N~ a 
with the fitting exponent a ~ 2.44, as shown in Fig(5Jb). 
With the emergence of a giant component of new phase, 
the tail of the distribution is elevated, but the size distri- 
bution for the remaining clusters still follows power-law. 
The underlying mechanism of such phenomenon is still 
an open question for us. 

To determine the critical size of the nucleus, A c , we 
compute the committor probability Pg, which is the 
probability of reaching the thermodynamic stable state 
before returning to the metastable state. The depen- 
dence of Pb on A is plotted in FigJ3Ja). As commonly 



FIG. 2: (Color online) (a) The average degree of nodes of new 
phase (Jz new ') as a function of A; (b) The size distribution of 
clusters of new phase, in which a power law distribution is 
present. Other parameters are N = 1000, (k) = 6, T = 2.59, 
and h = 0.7. 



reported in the literature [13, HH , the critical nucleus ap- 
pears at Pb(A c ) = 0.5, giving the critical nucleus size 
Xf FS — 474. The committor distribution at Af F s ex- 
hibits a peak at 0.5, of which 70% of spin configurations 
have Pb values within the range of 0.4 to 0.6 (see the 
inset of FigJHJa)), indicating that A is a proper order pa- 
rameter. 

Note that conventionally the nucleation threshold A c 
is usually estimated by using CNT [42TI441 ]. One can 
calculate the free energy change along the nucleation 
path, A-F(A), by using methods like umbrella sampling 
(US) Hi]. According to CNT, AF will bypass a max- 
imum at A = Ag , and the nucleation rate is given by 
v exp(— /3AF C ), where v is an attempt frequency. Here 
we have computed AF by using US, in which we have 
adopted a bias potential 0.1fcsT(A — A) 2 , with A being 
the center of each window. As shown in FigJ3Jb), the 
maximum in AF occurs at X^ s = 451, giving a free- 
energy barrier of AF C ~ 9lAksT. Clearly A^ gives 
a fairly good estimation of \ F FS . To calculate the nu- 
cleation rate, however, one has to obtain the attempt 
frequency v, which is not a trivial task. If we just 
set v = 1, we obtain a CNT prediction of a rate of 
2.02 x 10~ 40 MC step' 1 spin' 1 , which is 9 orders of magni- 
tude faster than that computed from FFS method. This 
level of disagreement in nucleation rate corresponds to an 
error in the free-energy barrier of about 24%. Since the 
accurate value of v is generally unavailable, we will use 
FFS method to calculate the nucleation rate throughout 
this paper. In addition, real nucleation pathway cannot 
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FIG. 3: (Color online) (a) The committor probability Pb as 
a function of A; The inset plots the committor distribution 
at Af FS . (b) The free energy AF as a function of A, in 
which the maximum in AF occurs at A^ s . (c) The logarithm 
of homogeneous nucleation rate InRnom (left axis), and the 
critical size of nucleus, A FFS and A^ s (right axis), obtained 
by FFS method and US method, respectively, as functions of 
h. Other parameters are the same as Fig(2] 



be obtained by conventional US method due to the use 
of a biased potential. 

We have also investigated how the nucleation rate and 
threshold depend on the external field h. In Figl3jc), 
InRnom, Af FS and are plotted as functions of h. 
The error-bars are obtained via 20 different network re- 
alizations and 10 independent FFS samplings. As ex- 
pected, lni?^f om increases monotonously with h, and 
Af and A^ both decrease with h. For large h, the 
difference between Af and X^ s becomes small. If h 
is large enough, one expects that the free-energy barrier 
will disappear, and nucleation will be not relevant. 



simulations to calculate Rnom and A c for different net- 
work size N. In particular, besides the BA-SFNs, we have 
also considered different network types, including SFNs 
with other scaling exponent 7 and homogeneous random 
networks (HoRNs) [46[. The networks are generated ac- 
cording to the Molloy-Reed model [47j : Each node is as- 
signed a random number of stubs k that are drawn from 
a specified degree distribution. This construction elimi- 
nates the degree correlations between neighboring nodes. 
We note here that the exponent 7 can be a measure of 
degree heterogeneity of the network, i.e., the smaller 7 
is, the degree distribution is more heterogeneous. In a 
HoRN, each node is equivalently connected to other (k) 
nodes, randomly selected from the whole network, such 
that no degree heterogeneity exists. By comparing the re- 
sults in SFNs with different 7 as well as that in HoRNs, 
one can on one hand, check the robustness of the system 
size effects, and on the other hand, investigate how the 
degree heterogeneity affects the nucleation process. 

Figure 2] shows the simulation results. All the param- 
eters are the same as in Fig. 2, expect that N varies from 
N = 500 to N — 3000. Interestingly, both \nR Hom 
and A c show very good linear dependences on the sys- 
tem size, i.e., lni?# om ~ — aN and A c ~ bN with a and 
b being positive constants. Obviously, in the thermody- 
namic limit N — > 00, we have Rnom and A c — > 00. 
This means that nucleation in these systems is not rele- 
vant in the thermodynamic limit, and only finite-size sys- 
tems are of interest. As shown in Fig. 4, for the network 
types considered here, qualitative behaviors are the same. 
Quantitatively, with increasing 7, the line slope becomes 
smaller, RHom becomes larger and A c gets smaller. Since 
larger 7 corresponds to more homogeneous degree dis- 
tribution, these results indicate that the degree hetero- 
geneity is unfavorable to nucleation. This is consistent 
with the nucleation pathway as shown in Fig. 2: In a het- 
erogenous network, those hub nodes are difficult to flip, 
making the nucleation difficult. 

In the following, we will show that the system-size ef- 
fects can be qualitatively understood by CNT and simple 
mean-field (MF) analysis. According to CNT, the forma- 
tion of a nucleus lies in two competing factors: the energy 
cost of creating a new up spin which favors the growth of 
the nucleus, and an opposing factor which is due to the 
creation of new interfaces between up and down spins. 
The change in the free energy may be written as [42j-[44| 



AF(A) = -2hX + aX, 



(3) 



B. Homogeneous Nucleation: System-size effects 

According to Fig [31 one finds that nearly half of the 
nodes must be inverted to achieve nucleation. This means 
that for a large network, nucleation is very difficult to 
occur. An interesting question thus arises: How the nu- 
cleation rate and threshold depend on the network size? 

To answer this question, we have performed extensive 



where a denotes the effective interfacial free energy, 
which may depend on T, h, and N. Since interfacial 
interactions arise from up spins inside the nucleus and 
down spins outside it, one may write a = 2JK out by 
neglecting entropy effects (zero-temperature approxima- 
tion), where K out is the average number of neighboring 
down-spin nodes that an up-spin node has. Using MF 
approximation, one has K out = (k) (1 — X/N). Inserting 
this relation to Eq[3]and maximizing AF with respect to 
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FIG. 4: (Color online) (a) The logarithm of homogeneous 
nucleation rate In Ruom and (b) the critical size of nucleus A c 
as functions of the network size N in SFNs with different 7 
and in HoRNs. Other parameters are the same as Fig[2] 



A, we have 



MF = J(k) h 



2J(k) 

and the free-energy barrier 

(J(k) -h) 2 N 



AF MF = 



2J(k) 



(4) 



(5) 



Clearly, both \^f F and AF^ IF linearly increase with N 
if other parameters are fixed. Therefore, the linear rela- 
tionships shown in FigH]are essentially analogous to the 
behavior of a mean-field network. Quantitatively, how- 
ever, the MF analysis fails to predict the line slopes in 
FiglU This can be understood because the approxima- 
tions are so crude, wherein important aspects such as 
network heterogeneity and entropy effects have not been 
accounted for. A rigid analysis is not a trivial task and 
beyond the scope of the present work. 



C. Heterogeneous nucleation 

In practice, most nucleation events that occur in na- 
ture are heterogeneous, i.e., impurities of the new phase 
are initially present. It is well known that impurities can 
increase nucleation rate by as much as several orders of 
magnitude. In our model, impurities are introduced by 



fixing some nodes in up-spin state. We are interested 
in how the number w of impurity nodes and the way of 
adding impurities would affect the nucleation rate. The 
first way of adding impurities we use is that impurity 
nodes are selected in a random fashion. Figure [5ja) gives 



the simulation results of In 



as a function of w in 



different network models, where Ruet are the rates of het- 
erogeneous nucleation. As expected, nucleation becomes 
faster in the presence of random impurities no matter 
which kind of network model is applied. It seems that in 
FiglSJa) all data collapse and exhibit a linear dependence 
on w, with the fitting slope 3.33. This means that each 
additional random impurity can lead to the increase of 
the rate by more than one order of magnitude. For the 
second way, we select w nodes with most highly degree as 
the impurity nodes, termed as targeted impurities. Strik- 
ingly, such a targeted scheme is much more effective in 
increasing the nucleation rate than random one, as shown 
in FigJSJb). For example, for SFNs with 7 = 3, one sin- 
gle targeted impurity can increase the rate by about 36 
orders of magnitude. 

As in Sec lIII Bl below we will also give a MF analy- 
sis of the heterogeneous nucleation, which qualitatively 
agrees with the simulation results. Each impurity node 
contributes an additional term to the free energy bar- 
rier, which can, under zero-temperature approximation, 
be written as the product of — 2 J and the expected 
degree of the impurity node. For random impurities, 
each impurity node has an expected degree (k), yielding 
the term —2J(k). Thus, the resulting free-energy bar- 
rier of the heterogenous nucleation becomes AiF^ et = 
AF H om _ 2 j ^ w> where AF H om is ^ f ree _ energy bar- 
rier of homogeneous nucleation. According to CNT, one 
obtains 



In 



Rffe 

R 



Ho 



2 J (k) 
kuT 



(6) 



Therefore, nucleation with random impurities is always 
faster than without impurity, and In {^ ^ ct J should vary 

linearly with w. The theoretical estimate of the slope is 
given by 2 J (k) /ksT = 4.63, approximately consistent 
with the simulation one. Given the simple nature of the 
above approximation the agreement is satisfactory. For 
the targeted way, a similar treatment to the former case 
can also be executed, except that (k) should be replaced 
by (k) w , where (k) w is the average degree of the w tar- 
geted nodes. After simple calculations, we can obtain 
(k) w = (k) (■^)^ rT - This leads to a free-energy barrier, 
- 2J(fc) N^w^, and 



A 2 F c Het = A 2 F^ om 



In 



( R 



Het 



\RHo 



2J{k) 1 

jv 1- 1 w t- 

k R T 



(7) 



Compared with EqJBl besides the presence of an addi- 

1 

tional size-dependent factor of N-/- 1 , the w-dependent 

T-2 

factor becomes w- 1 rather than w. For homogeneous 
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FIG. 5: (Color online) In ( ^ Hl!t ) as a function of the number 



\Rh 

of impurity nodes w. (a) corresponds to the case of random 
impurities, while (b) and (c) to the case of targeted impurities, 
(c) plots In f j^p^ J vs w in double logarithmic coordinates. 
All dotted lines are drawn by linear fitting. Other parameters 
are the same as FigO except for h = 0.2. 



network, i.e. 7 — > 00, one obtains N1- 1 — > 1 and 
7—2 

w — > w, and EqJT] is thus equivalent to EqjBJ In 

FiglSJc), we plot the simulation results of In f j^jp^ J as 

a function of w in double logarithmic coordinates, where 
linear dependences are apparent, in agreement with the 
result of Eq(?l The fitting slopes fi and intercepts K for 
7 = 2.55,3.0,5.0 in Figl^c) are fi = 0.49,0.62,0.93, 
and k = 2.25, 1.94, 1.08, respectively, while our ana- 
lytical estimates are \i = = 0.35, 0.50, 0.75, and 



= log 



2J{k)N-i- 
k B T 



= 2.60,2.16,1.42, respectively. 



Simulations and analysis give the same trends of the val- 
ues of /I, k and 7. 



rate is size-dependent, and it decreases exponentially 
with the network size, resulting in that nucleation only 
occurs at finite-size systems. As we already know, such 
system-size dependence does not exist in two-dimensional 
regular lattices [26]. We have also studied nucleation of 
Ising system in regular networks where each node is con- 
nected to its fc-nearest neighbors (here we only consider 
the case of sparse networks, i.e., k <C N), and found 
that the rate is almost independent of network size (re- 
sults not shown). Therefore, nucleation process in SFNs 
and HoRNs are quite different from that in regular lat- 
tices or networks. Such differences may originate from 
the infinite-dimensional properties of SFNs and HoRNs, 
wherein the average path distance is rather small, render- 
ing the system's behavior analogous to that of a mean- 
field network. An interesting situation arises when one 
considers Watts-Strogatz small-world network, which is 
constructed by randomly rewiring each link of a regular 
network with the probability p [48j . With the increment 
of p from to 1, the resulting network changes from a 
regular network to a completely random one. As men- 
tioned above, for the nucleation process, no system-size 
effects exist for p = 0, while system-size dependence ex- 
ists for p — 1. One may naturally ask: How does the 
crossover happens when p changes from to 1, and what 
is the physical significance of such a transition? This 
question surely deserves further investigations and may 
be the content of a future publication. 

In summary, we have studied homogeneous and hetero- 
geneous nucleation of Ising model in SFNs, by using FFS 
method. For homogeneous nucleation, we find that the 
formation of new phase starts from nodes with smaller 
degree, while nodes with higher degree are more stable. 
Extensive simulations show that the nucleation rate de- 
creases exponentially with the network size N, and the 
nucleation threshold increases linearly with N, indicat- 
ing that nucleation in these systems are not relevant in 
the thermodynamic limit. For heterogeneous nucleation, 
target impurities are shown to be much more efficient to 
enhance the nucleation rate than random ones. Simple 
MF analysis is also present to qualitatively illustrate the 
simulation results. Our study may provide valuable un- 
derstanding how first-order phase transition takes place 
in network-organized systems, and how to effectively con- 
trol the rate of such a process. 
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